Poincare recurrence and spectral cascades in 3D quantum turbulence 
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The time evolution of tfie ground state wave function of a zero-temperature Bose-Einstein con- 
densate (BEC) is well described by the Hamiltonian Gross-Pitaevskii (GP) equation. Using a set of 
appropriately interleaved unitary collision-stream operators, a quantum lattice gas algorithm is de- 
vised which on taking moments recovers the Gross-Pitaevskii (GP) equation in the diffusion ordering 
(time scales as length^). Unexpectedly, there is a class of initial conditions in which their Poincare 
recurrence is extremely short. As expected, the Poincare recurrence time scales with diffusion or- 
dering as the grid is increased. The spectral results of Yepez et. al. [T] for quantum turbulence are 
revised and it is found that it is the compressible kinetic energy spectrum that exhibits 3 strong 
cascade regions: a small-/c classical Kolmogorov k~^^^, a steep semi-classical cascade region, and a 
large- k quantum vortex spectrum k~^ . For most evolution times the incompressible kinetic energy 
spectrum exhibits a somewhat robust quantum vortex spectrum of for an extended range in 
k with a k~^'^^ spectrum for intermediate k. For linear vortices of winding number 1 there is an 
intermittent loss of the quantum vortex cascade with its signature seen in the time evolution of the 
kinetic energy Ekin{t), the loss of the quantum vortex spectrum k~^ spectrum in the incompress- 
ible kinetic energy spectrum as well as the minimilization of the vortex core isosurfaces that would 
totally inhibit the Kelvin wave vortex cascade. In the time intervals around these intermittencies 
the incompressible kinetic energy also exhibits a multi-cascade spectrum. 



I. INTRODUCTION 

The ground state of the many body wave function of 
a zero-temperature BEC is well described by the single 
particle wave function cp since all the bosons are in the 
same state. Upon appropriate normalization, the evolu- 
tion of this single particle wave function is governed by 
the Gross-Pitaevskii (GP) equation [3] 



ihdt(p = 



2m 



vV + (^I^I'-m)^, 



(la) 



for a spinless condensate where g is the nonlinear cou- 
pling representing the 5— wave scattering strength of the 
weak bosonic interactions in the mean-field approxima- 
tion and /i is the chemical potential. For numerical pur- 
poses it is convenient to rewrite (la) as 



idt(p 



-vV + «(^l^l'-i)^, 



(lb) 



where we consider a condensate comprised of particles of 
mass m = 1/2 in natural units where h = 1 and where 
a is introduced as a radial scaling parameter (for tuning 
the numerical resolution of the vortex cores). The "single 
particle" wave function cp is also the order parameter for 
the scalar BEC with \(p\ at each quantum vortex 
core (a topological singularity). 

The GP system is Hamiltonian. Now it is well known 
that for every Hamiltonian system there is a Poincare 
recurrence of the initial conditions, in that the dynam- 
ics will eventually return arbitrarily close to their ini- 
tial state — provided the dynamics are followed for a suf- 
ficiently long time [4 . For nearly all continuous Hamilto- 
nian systems, this Poincare time is, however, effectively 



infinite. Here we will show that there exists a special class 
of initial conditions for which this Poincare recurrence 
time is remarkably short and accessible with present day 
supercomputer resources. To define this class we follow 
Nore et. al., [5] and split the total conserved energy into 
several components 

^TOT = ^kin(t) + ^qu(t) + ^int(^) = COUSt. (2) 

where the kinetic energy £^kin(^), quantum energy £^qu(^) 
and internal energy £^int(^) are defined by 



Eintit) 



(3) 



on using the Madelung transformation 

(/p = V^e^^/^ with v = \/0 (4) 

which relates the condensate wave function cp to the BEC 
fluid density p and velocity v. We shall find that for 
initial conditions with 



^int(O) <^kin(0),^qu(0) 



(5) 



the corresponding Poincare recurrence time is extremely 
short. 

We have performed quantum turbulence simulations 
on the GP equation on grids of 5760^ [1 and re- 
ported on the existence of 3 distinct energy cascade re- 
gions. Due to an oversight in these earlier computations, 
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we attributed this triple cascade to the incompressible 
kinetic energy spectrum when, in fact, it should be basi- 
cally attributed to the total (and hence the compressible 
) kinetic energy spectrum. This is verified by simula- 
tions on lattice grids of 3072^ with details presented in 
the Appendix. In these new simulations one sees that 
the compressible kinetic energy spectrum exhibits 3 dis- 
tinct spectral regions: a small k classical Kolmogorov- 
like k~^^^ spectrum, a steep semi-classical spectrum 
with a non-universal exponent a > 6, followed by a ro- 
bust quantum energy spectral region with k~^ for large 
k. The incompressible kinetic is more complex. There is 
a very robust k~^ spectrum for a long wave number range 
extending to the maximally resolved k in the simulation. 
However, around the wave number k where there is a 
switch from the semi-classical to quantum vortex energy 
spectrum in the compressible spectrum, the incompress- 
ible kinetic energy spectrum does exhibit a well-defined 
k range in which the spectrum is k~^^^^. 

The identification of the k~^ and k~^^^^ spectra is 
strongly debated. Since we have many quantum vor- 
tices present we will call the large- /c spectrum the quan- 
tum vortex spectrum. As Nore et. al. [5 have pointed 
out, the energy spectrum of a single stationary linear 
quantum vortex is readily determined analytically and 
it exhibits only one long single spectrum of all the 
way to the maximum k. It is readily shown that this 
holds for both a 2D and a 3D single vortex. For either 
a Taylor- Green or Berloff vortex, all the kinetic energy 
is incompressible - there is no compressible energy spec- 
trum. Of course, a k~^ incompressible kinetic energy 
spectrum that is found in the dynamical evolution of the 
GP wave function, does not necessarily imply a simple 
quantum vortex. The k~^^^^ incompressible kinetic en- 
ergy spectrum is also quite interesting. There is also a 
strong debate on the spectrum of quantum Kelvin waves 
on a quantum vortex. Kozik & Svistunov [9 argue for 
a k~^^l^ spectrum in the wave-action spectrum. It has 
been shown [10] that in some regimes, the wave-action 
spectrum will have the same exponent as the kinetic en- 
ergy spectra, and so it is not impossible that the k~^^^^ 
spectrum we observe could be associated with the quan- 
tum Kelvin wave cascade. As regards the full kinetic en- 
ergy spectrum, one can identify a classical Kolmogorov 
^-5/3 spectrum for small k. It is interesting to note that 
in subgrid scale turbulence closure modeling of sonic jet 
flows in compressible classical turbulence, it is deemed 
critical [11] that the closure model reproduce the Kol- 
mogorov scaling for the total kinetic energy. 

In most of the simulations considered here we run on 
grids of 1200^ since the similar physics (Poincare recur- 
rence and the associated phenomena) on a 3072^ would 
require 6.6 times more iterations. What we shall find is a 
robust k~^ incompressible kinetic energy which is inter- 
mittently destroyed together with a minimization in the 
vortex core isosurfaces. Away from these intermitten- 
cies one sees a dual kinetic energy spectrum. The total 
kinetic energy spectrum exhibits three different spectral 



regions - although somewhat weaker than the three spec- 
tral regions in the 3072^ simulations. 

In Section 2 we describe our unitary quantum lattice 
gas (QLG) algorithm which is extremely well parallelized 
on supercomputers (no performance saturation seen with 
scaling runs tested to 216000 cores )and which will run on 
quantum computers when available. Moreover QLG also 
has a low memory footprint. In Section 3 we discuss the 
Poincare recurrence for a set of initial conditions satisfy- 
ing Eq. (5). The signature of this recurrence is readily 
apparent on plotting the time evolution of the kinetic 
energy £^kin(^) and quantum energy E^u{t). Initially the 
rectilinear vortices are very well localized. As a result, 
for nearly all regions in coordinate space the BEG density 
p is constant. From Eq. (3) this translates into very low 
initial quantum kinetic energy £^qu(0) <C £^kin(0) with 
^comp(O) <C £^incomp(0) • The Poiucarc recurrence is fur- 
ther established by looking at the vortex core isosurfaces. 
These simulations are run on grids of 1200^. In the con- 
tinuum limit that yields Eq. (1), our QLG algorithm 
follows diffusion ordering. This has been verified by in- 
creasing the grid from to L2 and seeing the Poincare 
recurrence time increase from Tpoin to {L2/ Li)'^Tpoin 
(grids used were 512^, 960^, 1024^ and 1200^). In Section 
4 we shall discuss the somewhat unexpected intermittent 
loss of the quantum vortex spectrum in these Poincare 
recurrence simulations. We find on examining the incom- 
pressible kinetic energy spectrum that the pronounced 
k~^ spectrum is destroyed for approximately 2000 time 
steps. A precursor for the loss of this vortex cascade 
is identified in the time evolution of the kinetic energy 
^kin(^)- In Section 4 we also see how the topological 
changes in the quantum vortex cores (during this inter- 
mittent loss of energy spectrum) are also registered in 
the minimization of the vortex core isosurfaces. Clearly, 
this minimization will preclude a quantum Kelvin wave 
cascade from occurring along the vortex core. It is inter- 
esting to note [12 that in 2D GP quantum simulations, 
where the quantum vortex cores are just point singulari- 
ties, one also sees the k~^ large-k incompressible kinetic 
energy spectrum provided there were point vortex cores 
in the dynamics. When the point vortices are annihilated 
during the GP dynamics, this k~^ spectrum is lost. 

We here briefly put our QLG method into the context 
of the mesoscopic literature. In computer simulations of 
nonlinear physics, it is desirable to develop algorithms 
that scale strongly with increasing number of processors. 
Unfortunately, nearly all direct solution methods do not 
scale well because of the non-local nature of the equa- 
tions. Lattice Boltzmann and lattice gas algorithms [13] 
move to a mesoscopic level by either solving a simpli- 
fied discretized collisional Boltzmann equation or by a 
simplified particle collide-stream set of rules on a lattice. 
Constraints are imposed so that in the macroscopic limit 
one recovers the original system of equations. What is 
gained is a very local representation that is ideally par- 
allelized. One of the earliest quantum lattice Boltzmann 
schemes was by Succi and Benzi [14|. Succi [15^ then 
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extended his complex distribution function approach to 
examine the GP equation and extended the algorithm 
to 2D and 3D [16l [17] . A quantum lattice gas approach 
based on qubits, similar to what is proposed here, was 
initiated by Yepez [18^ J. 9 and extended by Yepez and 
Boghosian [20]. Quantum lattice algorithms were also 
examined by Meyer pT. As our qubit lattice represen- 
tation is a mesoscopic representation of the GP equation 
it is critical that it be benchmarked against exact solu- 
tions. In particular, our codes have been benchmarked 
against the ID soliton-soliton collisions of NLS in which 
the solitons retain their exact shape and speed together 
with a phase-induced shift due to the collision [22 . We 
have also benchmarked our codes for ID vector soliton- 
soliton collisions in which one can have inelastic collisions 
for very specifically chosen initial amplitudes [23^ ^4 . Fi- 
nally we comment on the Poincare recurrence time of ID 
NLS. At first glance it would seem that the Poincare re- 
currence time for such a continuum Hamiltonian system 
would be effectively infinite. However there were some 
hints in Tracy et. al. [2F that the solution manifold 
for ID NLS has an underlying finite-dimensional struc- 
ture. Thus, under certain initial conditions, one could 
find short Poincare recurrence times for ID NLS. What 
is more unexpected is that we have found a class of ini- 
tial conditions which yield short Poincare recurrence for 
3D quantum turbulence. Outside this class, the Poincare 
recurrence time is essentially infinite. 



II. QUANTUM LATTICE GAS ALGORITHM 
FOR THE GP EQUATION 

To recover the scalar GP equation in the continuum 
limit, we consider two qubits at each lattice site. Of 
the four possible states |00), |01), |10) and |11), we need 
consider just the complex amplitudes a and {3 for the 
states 1 01), 1 10) respectively. At each position x of the 
cubic lattice we introduce the two-spinor 



a(ic, t) 
/3ix,t) 



(6) 



and construct the evolution operator U[Q] - consisting of 
an appropriate sequence of non-commuting unitary colli- 
sion and streaming operators - so that in the continuum 
limit the two spinor equation 



+ At) = u[n]iij{x,t). 



(7) 



will reduce to the GP equation for the 1-particle boson 
wave function (p under the projection 



(1,1) •^ = ¥' 



(8) 



The unitary collision operator C that entangles locally 
the amplitudes a and /3 is chosen to be the square-root-of 
swap: 



(J — g«f crcc(l-cra;) 



1-i 1+i 

2 2 

1+i 1-i j 7 

2 2 



(9) 



where the a are the Pauli spin matrices 



on _ (0 -i 

10 ~ U 



is just the swap gate since 
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(10) 



(11) 



The streaming operators shift just one of these am- 
plitudes at X to a neighboring lattice point at x + Ax: 



£C,0 



a{x^ t) 
p{x,t) 

a{x^ t) 
p{x,t) 



a{x + Ax, t) 
p{x,t) 

a{x^ t) 
f3{x + Ax,t) 



(12a) 
(12b) 



The subscript 7 = on the streaming operator 5'a£c,7 
refers to shifting the amplitude a while the subscript 7 = 
1 refers to shifting the amplitude (3. In terms of the Paul 
spin matrices, the streaming operators can be written in 
the form 



Axe 



Sax,i = n + e 



Ax6 



(13) 



where n = (1 — crz)/2^ n = (1 + az)/2. It should be 
noted that the collision and streaming operators do not 
commute: [C, S] ^ 0. 

We now consider the following interleaved sequence of 
unitary collision and streaming operators 

Jx^ — S—Ax,^C Sax,^C (14) 

Since \Ax\ <C 1 and = I, J^^ = I ^ 0{Ax), where 
/ is the identity operator. We first consider the effect of 
the evolution operator Uj[Q{x)] 



u^[n{x) 



t2 t2 t2 -i£^n{x) 
x'^^ y^"^ z-^^ "> 



(15) 



acting on the 7 component of the 2-spinor ?/^. Here £ <C 1 
is a perturbative parameter and O is a function to be 
specified later. 

Using perturbation theory, it can be shown that the 
time advancement of 



V^(ic,t + At) = U^[^]\lj{x,t). 

yields 

V^(ic, t + At) = iIj{x, t) - ie^ 
(-l)V 



(16) 



1 



':(Jx^' 



tlj{x,t)- 



-(a^+a,)VXa^,t) + (9(s^), 



(17) 



with 7 = or 1 and Aa; = 0(e). Since the order term 
in (17) changes sign with 7, one can eliminate this term 



by introducing the symmetrized evolution operator 



Up] = Ui 



Uo 



(18) 
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rather than just U^. 

Under diffusion ordering, At = 0{e'^) and Ax = 0{e)^ 
the evolution equation 

t^At) = U[n{x)] iIj{x, t) (19) 

leads to a representation of the spinor equation 

1 



idtil^{x,t) 



^{x,t)^0{e'), 



(20) 



where the function Vt is still arbitrary. To recover the 
scalar GP equation, one simply rescales the spatial grid 
V — > a~^V, contracts the 2-component field to the 
(scalar) BEG wave function (f 

(/p = (l,l)-^ = a + /3 (21) 

and chooses Vt = g\(p\'^ — 1 : 

idt^ = -VV + cl{9\^\^ - 1)^ + 0{e^). (22) 

One of the beauties of the QLG algorithm, Eqs. (20), 
(19), (16) and (15), is its ideal parallelization on cur- 
rent day supercomputers (as well its direct application 
to quantum computers - once available - due to its uni- 
tary representation). Indeed because the collision oper- 
ator is purely local and the streaming operator requires 
shifting data just to nearest neighbors, we have seen no 
saturation of the parallelization to over 163 000 cores 
on BlueGene / Intrepid (Argonne) and on over 216 000 
cores on Jaguarpf (Gray XT-5 at Oak Ridge National 
Laboratory) — these being the maximum number of cores 
available to us currently. 



III. POINCARE RECURRENCE FOR A CLASS 
OF INITIAL CONDITIONS 

We consider straight line vortices as the initial condi- 
tions for our simulations of the GP Eq. (22). Unlike a 
classical vortex, a quantum vortex is a topological singu- 
larity with the wave function (/:? = at the core singular- 
ity while the velocity Eq. (4), diverges there. Using 
Pade approximants [26 , on the GP Eq. (22) one can de- 
termine an asymptotic steady state straight line vortex. 
For winding number n = 1, using cylindrical polar coor- 
dinates (r, (/), z), such a vortex that lies along the 2;-axis 
(and centered at the origin) is given by 



I Ilar2(12 + ar2) 
^/[384 + ar2(128 + lla^ 



(23) 



with \(p\ 1/v^ as r ^ oo, and \(p\ ~ r^/a/g as r ^ 0. 

The coherence length ^, for this single line vortex, is 
typically defined as the distance from the core singularity 
to the position at which the absolute value of the wave 
function approaches its asymptotic value 



Quantitatively, the coherence length is defined from the 
solution of an isolated line vortex for a boundary value 
problem of the GP equation. Here, we study quantum 
turbulence with many vortices interacting under periodic 
boundary conditions. Alternatively, following Nore et. 
al. [5 , the coherence length can be defined from a linear 
perturbation dispersion relation about a uniform density. 
Qualitatively, in quantum turbulence with its many in- 
teracting quantum vortices, it is tempting to define a co- 
herence length by replacing the background asymptotic 
density po by the spatially averaged mean BEG density 
< Po >5 However, in our simulations we rescale 



the initial wave function, Eq.(23) so that our initial con- 



dition are far from a quasi-steady state solution to the 
GP Eq.(22) and so leads to turbulence more rapidly, es- 



-1/2 



(24) 



pecially on very large grids. Under these conditions we 
feel that the standard idea of coherence length does not 
readily apply. 

For winding number n = 1 there is a phase change 
of 27T in a closed circuit about the core singularity. To 
enforce periodic boundary conditions we consider a set of 
four vortices parallel to the z-axis, with a similar set of 
four vortices with axes parallel to the ?/-axis and x-axis. 
The total initial wave function is a product [4] of these 
twelve straight vortices and the initial quantum vortex 
core singularities are shown in Fig. [l] along with the 
phases on the boundaries of the lattice. One immediately 
sees that at those points where the vortex isosurfaces 
intersect the plane walls there is a phase change from 
(j) = (color blue) to = 27r (color red). Clearly there 
are 12 such points. Moreover, the sense of the BEG fiuid 
rotation about the core singularity can be gleaned from 
the rotation sense of ^ = to = 27r. Also shown in Fig. 
[T]is a second set of initial conditions considered: a set of 
12 doubly degenerate line vortices (essentially winding 
number of n = 2 with \(p\ ^ r'^a/g as r ^ 0) with phase 
change of 47r around each core singularity. The double 
degeneracy is evident from the phase information on the 
walls about these core singularity intersections and that 
this degeneracy is easily broken in the time evolution 
of the BEG as it is simply a confiuence of two winding 
number n = 1 singularities with energies significantly 
lower than that of the degenerate state. 

We plot isosurfaces very close to the vortex core sin- 
gularity. The parameters and rescaling of the wave func- 
tion are so chosen that Eint{0) <C Ekin{^)^ Equ{0)^ with 
a{0) = /3{0) = (p{0)/2. The early time evolution of the 
vortex core singularities are shown in Fig. [2] (for winding 
number n = 1) and in Fig. [3](for winding number n = 2). 
Here the phase information both on the walls (with the 
green core singularity) and on the singularity core (with 
grey wall background) are plotted. The 27r phase change 
is evident when the singularity core intersects the walls 
[Figs. [2|a), (c) and (e)] while on close inspection one also 
sees the 27r phase change along the singularity core itself. 
Fig. ^h), (d) and (f). Kevlin waves along the cores are 
clearly seen as are reconnect ions and vortex loops. 

As expected the vortex core singularity structures for 
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(a) winding no. n=l 




(b) winding no. n=2 



FIG. 1. The initial isosurfaces of the 12 vortex core singularities 
for (a) winding number n = 1, and (b) winding number n = 2. 
(a) For winding number n = 1, the core isosurface shown here is 
for \ip\ = 0.13|(/?|max since for lower isosurface values these initial 
cores would be too faint, (b) For winding number n = 2, the core 
singularities are shown for isosurface \<^\ = 0.008|(/:?|max The initial 
location of the 12 vortex cores for winding numbers 1 and 2 are 
different. On the walls are shown the corresponding phases cf), with 
(f) = in blue and (j) = 2tt \n red. The core singularity intercept 
with the walls acts like a branch point in the phase plane, with a 
branch cut joining the branch points. Grid 1200^. 



winding number n = 2 (Fig. |3| are much more pro- 
nounced than for those core structures for winding num- 
bern = 1 (Fig. [2|. 

The time evolution of the total energies, Eq. (3), is 



shown in Fig. |4]for vortex cores with (a) winding number 
n = 1 and (b) with winding number n = 2. Note that the 
initial quantum vortices are so highly localized in space 
that the BEG density p ~ const throughout the cubic 
lattice. Thus V s/p ~ throughout the lattice so that 

Equ{^) <C Ekin{^) and Ecomp{^) ^ Eincomp{0). 

For quantum vortex cores with winding number n = 1, 
from the time evolution of the kinetic and quantum ener- 
gies one sees that £^kin(0) ~ £^kin(^ = 115000) ^ £^kin(^ = 
230000) ^ • • • . The situation is a bit more complex for 
cores with winding number n = 2 in that the subsequent 
peaking in £^kin at t = 115000 and t = 230000 is some- 
what reduced from its initial value. This return of £^kin(^) 
and £^qu(^) close to their initial values is suggestive of 
a Poincare-recurrence, particularly as the GP equation 
is a Hamiltonian system. To examine this we consider 
the corresponding vortex core isosurfaces for the differ- 
ent winding numbers at times t = 115000, Fig. [5] . As 
expected, for winding number n = 2, the confluence of 
the two singularities is broken with a slight split into two 
adjacent non-degenerate winding number n = 1 singular- 
ities. The vortex core singularity isosurfaces for winding 
number 1 at t = 230000 are shown in Fig. |6] where it 
is clear that the Poincare recurrence time on a 1200^ 
grid is Tpoin = 230000. On comparing Fig. [ijwith Fig. 
|6] one does see some slight variations in the wall phase 
information from that at t = as well as some small 
scale isosurface loop and vortex features now appearing 
for winding number 2. 

It is of some interest to note that the point inversion 
of the initial conditions at Tpoin/'^ for the BEG wave 
function is similar to what can be seen in the discrete 
Arnold cat map for some^ but not all, pixel resolutions 
of the initial state - see Fig. [7] 

The robustness of the Poincare recurrence is further ex- 
hibited on simulating the evolution of 48 quantum vortex 
core singularities, Fig. [Sj 

Further, we find (see Table I) that the Poincare re- 
currence time, Tpoin^ scales with diffusion ordering: i.e., 
as the lattice grid scales from Li L2 the Poincare 
time rescales from Tpoin (^2/^i)^^Poin- It should 
be stressed that there is nothing per se in our quantum 
algorithm, Eq. (19), that enforces diffusion ordering at 
the GP level - to recover the GP equation it is critical 
to have diffusion ordering in our qubit unitary algorithm 
and this is achieved by specific choices of parameters in 
the algorithm itself. Thus the Poincare recurrence fol- 
lows physics ordering rather than the naively expected 
lattice ordering (Ls/Li)^ for 3D. 



IV. INTERMITTENT LOSS OF THE LARGE-K 
INCOMPRESSIBLE KINETIC ENERGY 
SPECTRAL CASCADE 

An unexpected feature of our simulations is the inter- 
mittent loss of the incompressible kinetic energy spectral 
cascade for simple sets of initial line vortices with wind- 
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(e) n=l8itt = 99000 (f ) n = 1 at t = 99000 phase colored 



FIG. 2. The evolution of the isosurface core singularity (defined by \(^\ = 0.008|(y:?|max) for winding number n = 1 at times (a) t = 6000, 
(b) also at t = 6000 but zoomed in perspective (same orientation) and now with phase information on the core itself; it is evident that 
there is a 27r phase change around each vortex with vortex bending, reconnection, ring and blob formations. Similarly (c) and (d) for time 
t = 25000 and (e) and (f) for time t = 99000. Phase coding : cf) = in blue, cf) = 27t in red. Grid 1200^. 
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FIG. 3. The corresponding isosurfaces of the core singularities for winding number n = 2 at t = 6000 (with zoomed-in perspective for 
phase-coded core singularity), t = 25000 and t = 99000. Phase coding : (j) = in blue, = 27r in red. Grid 1200"^ 
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FIG. 4. The time evolution of the £^kin(^) (in blue) and Ec^u(t) 
(in red) for < f < 250000 for (a) Winding Number n = 1, and (b) 
Winding Number n = 2. Throughout the run the unitary algorithm 
very well preserves the total energy conservation £^tot = const. 
The internal energy is negligible. Note that £^kin(0) ~ ^kin(^ = 
115000) ^ E^i^{t = 230000) ^ ••• for vortex cores with winding 
number n = 1. For vortex cores with winding number n = 2, 
Ehin(t = 115000) and to a somewhat lesser extent £^kin(^ = 230000) 
tend to £^kin(^ = 0)- 




(a) n=l,t = 115000 




(b) n 



115000 



TABLE L The scaling of the Poincare recurrence time, Tpoin 
with the lattice grid L assuming diffusion scaling and comparing 
this theoretical scaling with simulation results. The scalings are 
normalized to a lattice grid of 512^. 



Grid 


Diffusion Scaling (theory) 


Simulation 


512^ 


41775 


41775 


960^ 


146850 


147000 


1024^ 


167100 


167600 


1200^ 


229477 


230000 



FIG. 5. The isosurface cores, with phase information on the 
walls, around the semi-Poincare period t ~ 115000 for (a) winding 
number 1, and (b) winding number 2. One sees a point inversion 
of the initial isosurface vortex coers. For visibility, the vortex core 
isosurface for winding number 1 is shown at |(^| = 0.064|(/?|max, 
while those for winding number 2 at \lp\ = 0.008|(y9|max- Phase 
coding : = in blue, </) = 27r in red. Grid 1200"^ 



ible and incompressible components 



ing number n = 1. This intermittency loss was first es- 
tablished from the time evolution of the incompressible 
kinetic energy spectrum and then verified by the topo- 
logical behavior of the quantum vortex core singularities 
during this time interval. To determine the incompress- 
ible kinetic energy spectra one first must decompose the 
density weighted velocity into its orthogonal compress- 



(25) 



with 



^ ' [^/p'^jincomp — 



(26) 
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(a) n = l,t = 230000 




(b) n = 2,t = 230000 



FIG. 6. The isosurface cores around the Poincare period t ~ 
230000 for (a) winding number 1, and (b) winding number 2. The 
confluence degeneracy of winding number 2 singularities are now 
broken globally and there is considerable vortex loop generation 
locally. Phase coding : (f) = in blue, cf) = 27r in red. 



and thus obtain [3^ the compressible and incompressible 
kinetic energy spectra 

i?fe,„(t)=4™w+^fe?rw (27) 

Jo Jo 

For the initial 12 straight line quantum vortices, over 98% 
of the initial kinetic energy is incompressible, Fig. [9ja) 

By t > 8000 the kinetic energy spectra have rapidly 
reached their quasi-steady state forms. As can be seen 
in Figs, and (c), the compressible kinetic energy 




(a) initial image (b) iteration 80 




(c) point inversion symmetry, (d) Poincare recurrence at 
iteration 157 = Tp/2 iteration 314 = Tp 



FIG. 7. The existence of point inversion symmetry for a given 
pixel resolution of an image under the Arnold Cat map (x, y) — 
(2x -\- y,x -\- y)modl. This 2D map is invertible, area preserving, 
ergodic and mixing. If the resolution of the skyrmion, (a), is 313 x 
313, then at iteration 80 one has the image (b). At Tp/2 = 157, (c), 
one has the point inversion of the original skyrmion while at Tp = 

314 one has Poincare recurrence. However, for pixel resolution 

315 X 315 there is no point inversion symmetry and Tp = 120. (A 
skyrmion is a two-component topological defect with a vortex ring 
threaded by a quantum vortex line that then closes on the outer 
surface like an apple core). 



spectrum exhibits 3 cascade regions while the incom- 
pressible kinetic energy spectrum exhibits basically two 
slopes. It should be pointed out that we had an error in 
our FFTs in [1 and the spectra exhibited is basically the 
compressible spectrum — see also [27j. To examine the 
Poincare recurrence of initial conditions and the inter- 
mittent loss of the vortex spectrum, we consider grids of 
1200^ - for example, for grids of 3072^, the first Poincare 
recurrence would occur around t = 720 000 and require 
much wallclock time and cores with little new physics 
insights. For t > 8000 the compressible kinetic energy 
spectrum has 3 distinct cascade regions: a classi- 
cal Kolmogorov-like energy cascade for 15 < k < 35 
with exponent a ~ 5/3, a steep semi-classical cascade 
for 40 < /c < 70 and a quite long quantum vortex spec- 
trum for 100 < k < 450 with exponent a ~ 3. On the 
otherhand, the incompressible kinetic energy spectrum 
exhibits a dual cascade over nearly all wave number range 
with exponents 7 ^ 3.7 and 7 = 3.0. Since for t > 8000, 
> 3.5E~^{t) and from Fig. gone can read- 
ily see that the total kinetic energy spectrum will also 
exhibit these 3 cascades will little change to the corre- 
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(e) n = 2, t = , 48 linear vortices (f) n = 2, t = 230000 



FIG. 8. Evolution of quantum core singularities from an initial set of 48 straight line vortices with phase information on the walls, (a) 
winding number n = 1 at t = 0, (b) winding number n = 1 a^t t = 84000. (c) winding number n = 1 at t = 230000 = Tp showing only- 
small perturbative changes to the initial state given in (a), (d) winding number n = 1 a,t t = 115000 = Tp/2. The 27r phase changes at 
the core singularity intersections at the walls is very evident, (e) winding number n = 2 at t = 0. (f) The corresponding isosurface cores 
at t = 230000 = Tp for winding number n = 2. The wall phase is a simple perturbative change from that at t = 0, (e) - but much small 
scale vortex structures due to the confluent degeneracy . Phase coding : (f) = - blue, (j) = 27t - red. Grid 1200^ 
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500 1000 




500 1000 



500 1000 



(c) n=l ,t = 75000 

FIG. 9. Time snapshots of the incompressible (blue) and com- 
pressible (red) kinetic energy spectrum at times (a) t = 0, (b) t = 
20 000, and (c) t = 75 000. Initially nearly 99 % of the kinetic en- 
ergy is incompressible, while at t = 20 000 over 72% of the kinetic 
energy is now compressible. The spectra are quasi-steady state 
for 8000 < t < 81000. The compressible kinetic energy spectrum, 
k~^, exhibits 3 clear cascade regions: a classical Kolmogorov-like 
spectrum a ^ 5/3 for 15 < /c < 35, a semi-classical cascade for 
40 < /c < 70 and a quantum vortex spectrum with o; ~ 3 for 
100 < k < 450. The incompressible kinetic energy, on the other 
hand, exhibits a dual cascade spectrum k~^: for small k,a 3.7 
while for large k,a = 3.0. Grid 1200^ 



spending exponents. 

For 81400 < t < 84300, the kinetic energy spectra 
intermittently take the forms shown in Fig. 



10 The 



incompressible spectrum now takes on a quasi 3 cascade 
spectrum with a brief loss in the spectral region. 

If one examines the vortex core isosurfaces, Fig. [TTJ 
one sees that the topological structure of the vortices are 
minimal and thus do not permit vortex- wave interactions 



along these 'point-like' cores. 

It is interesting to note that the loss of the quantum 
vortex spectrum, as seen in the spectral change in 
the incompressible energy, can also be discerned from the 
time evolution of the kinetic and quantum energies. Fig. 
|4] Indeed, in Fig. [T2j we plot the detailed time evolution 
of Fig. [4] around the loss of vortex cascades (81400 < t < 
84300) and on the onset of the half-Poincare recurrence 
time Tp/2 = 115000. As is also evident from Fig. |4]and 
from simulations results, this loss of the spectrum 
reoccurs in the time interval 196400 < t < 199300 with 
the full Poincare recurrence of the initial conditions at 
Tp = 230000. However only the Poincare recurrence time 
carries over to the case of initial line vortices with winding 
number n = 2 - the loss of the vortex cascade is not 
seen as expected from the split in the degeneracy of the 
winding number 2 vortices (see Fig. [sb and Fig. Hb) and 
the time evolution of Ekin{t) and Eq^t) around (81400 < 
t < 84300). The vortex core isosurfaces are shown in Fig. 
[12] at t = 1000 to illustrate the degeneracy splitting of the 
winding number n = 2 vortices and at time t = 82000 
where some minimal vortices can be seen among those 
which can support Kelvin wave cascades. This strong 
minimalization of the vortex core isosurfaces at the loss 
of the vortex cascade is also seen at other grid resolution 
runs, with the time of occurrence scaling with diffusion 
ordering. However, it is not that surprising that this 
loss of the vortex cascades is not seen for higher winding 
numbers or for more complex initial conditions (e.g., like 
the 48- vortex case which still retains considerable vortex 
cores and loops at t = 84000) 

One possible explanation for these spectral results is 
based on the role of the spectrum. If one consid- 
ers an isolated quantum line vortex, its kinetic energy is 
all incompressible and this incompressible kinetic energy 
spectrum is, for all wave numbers, k~^ [5 . When we 
turn to the dynamics of quantum vortices we have found 
this k~^ spectrum for incompressible and compressible 
kinetic energy spectra as well as for the quantum energy 
spectrum for very large k. This k~^ spectrum for the 
incompressible kinetic energy spectrum is intermittently 
lost when vortex loops are no longer present in the turbu- 
lence, as seen in Fig. [lO] For more complex initial con- 
ditions, like winding number 2 vortices, the isosurfaces 
in the time interval 82200 < t < 84000 will retain vor- 
tex loops and the incompressible kinetic energy spectrum 
retains its k~^ spectrum. A very similar result has also 
been found in 2D GP turbulence simulations, where un- 
der initial random phase conditions and constant density 
(i.e., when there are initially no [point] vortices), many 
vortices are rapidly born. One again encounters both 
very short Poincare recurrence time Tp of that initial 
random phase as well as point inversion of these phases 
at Tp/2 - and during these times of no point vortices the 
k~^ incompressible kinetic energy spectrum is lost. 

There is a wave number k^ at which there is a spec- 
tral break in the incompressible kinetic energy spectrum, 
where for k < k^ has a spectral exponent a ~ 3.7 while 
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81400 82200 83 400 




(d) n = l ,t = 84000 (e) n=l ,t = 84200 (f ) n = l ,t = 85000 



FIG. 10. 6 snapshots of the incompressible (blue) and compressible (red) kinetic energy spectrum at times (a) t = 81 400, (b) t = 82 200, 
(c) t = 83 400, (d) t = 84 000, (e) t = 84 200, and (f) t = 85 000. In (b)-(d) there is a sharp drop in the incompressible energy spectrum 
for wave numbers k > 100, except for a very brief transient recovery around t ~ 83000. There is also a sharp cutoff in the compressible 
spectrum for k > 500. This signals the loss of the quantum vortex spectrum in this time interval. Around these intermittencies, the 
incompressible kinetic energy spectrum also exhibits a triple cascade k~^ with a ~ 3.7 for small fc, an o; ~ 6 for the intermediate cascade, 
and a ~ 3.0 for the large-fc quantum vortex spectrum. During the intermittency, the large-fc exponent increases to a noisy a ~ 5.2 as well 
as a steeped semi-classical intermediate exponent. Grid 1200^ 



for < k the exponent is a = 3. Moreover, at this 
wave number /c^, we also see the strong break in the com- 
pressible and quantum energy spectrum from the semi- 
classical regime with strong a > 6 to the ubiquitous 
a = 3. It is thus tempting to associate the wave number 
spectrum k^ < k with that of an isolated vortex - together 
with sound/shock waves with sharp density variations so 
that the compressible and quantum energy spectra have 
a k~^. The sharper spectra for all the energies in the re- 
gion k < k^ could possibly be attributed to either kelvin 
waves on the quantized vortices or a Saffman-like k~^ 
spectrum due to vorticity discontinuities. On grids 1200^ 
the coherence wave number k^ ~ 70, while on grids 3072^ 
the k^ ^ 300, from the simulations presented in the Ap- 
pendix. 



It is also very tempting to associate the small k- 
region in the total energy spectrum with the classical 
Kolmogorov-like spectrum. It is important to note that 
this k~^^^ spectrum is not seen in the incompressible ki- 
netic energy spectrum - but in the total (and basically the 
compressible) kinetic energy spectrum. It is interesting 
to see that in classical compressible turbulence studies 
[TT] , it is the full subgrid energy spectrum that is required 
to model the Kolmogorov k~^^^ - not the incompressible 
part of this subgrid energy spectrum. 



V. CONCLUSION 

A novel quantum unitary lattice gas algorithm is de- 
vised to solve the time evolution of the ground state wave 
function of a zero temperature BEG as given by the GP 
equation. We introduce 2 qubits per spatial node and 
concentrate on the 1 body sector. A particular inter- 
leaved sequence of unitary square root of swap and uni- 
tary streaming operators act on a 2-spinor state. Pa- 
rameters are so chosen that under diffusion ordering the 
zeroth moment of the 2-spinor state reduces, in the con- 
tinuum limit, to the scalar wave function given by the GP 
equation. We find a particular set of initial conditions for 
which the Poincare recurrence is surprisingly short for 
this Hamiltonian system. Since the Poincare recurrence 
time scales with diffusion ordering, as verified by our sim- 
ulations on the 2-spinor state, we have run only to grids of 
1200^. It is also seen that the compressible kinetic energy 
spectrum exhibits 3 distinct power law cascades with the 
small k cascade corresponding to the Kolmogorov cascade 
of classical fiuid turbulence while the large k spectrum 
corresponds to the vortex cores themselves. The incom- 
pressible kinetic energy spectrum exhibits predominantly 
a dual cascade with spectra k~^-^^ for small k and k~^-^ 
for large k. It is interesting to note that recently Kerr [28] 
has also reported on a k~^ spectrum: in particular, Kerr 
considered the reconnection of two antiparallel quantum 
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(c) n = l ,t = 82000 (d) n = l ,t = 88000 



FIG. 11. 4 snapshots of the vortex core singularity isosurfaces at times (a) t = 78 000, (b) t = 81 000, (c) t = 82 000, and (d) t = 88 
000. Initial conditions of 12 line vortices. The phase information (blue is (j) = 0, red is (j) = 27t on the vortex core singularities clearly 
shows the 27r phase change in circumnavigating the vortex core. Grid 1200^ 



vortices. Following the first reconnection, vortex waves 
propagate along the quantum vortex and increase in am- 
plitude so that secondary reconnection occurs with the 
formation of vortex rings. At this stage Kerr [28] finds 
an energy spectrum of with the interaction energy 
concentrated around the vortex cores where p ~ 0. We 
introduce this work to stress that a spectrum may 
not necessarily imply the spectrum from a simple isolated 
vortex which indeed does also have a spectrum. In 



our simulations, the intermittent loss of the vortex cas- 
cade is readily seen in the loss of the spectrum of the 
incompressible kinetic energy spectrum as well as in the 
topological minimilization of the vortex cores. 

This work was partially supported by the Air Force 
Office of Scientific Research and the Department of En- 
ergy. Computations were predominantly performed on 
the SGI Altix ICE at the ARL and ERDC DoD High 
Performance Computing Centers, with some computa- 
tions also performed at DoE NERSC facility. 
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(a) n = l , spectra for 78000 < t < 86000 (b) n=l , spectra for 106000 < t < 120000 



FIG. 12. Segments in the time evolution of the total kinetic energy E]^in{t) (blue) and quantum energy Equ(t) (red) around (a) the loss 
of vortex cascade around 81400 <t < 84300, and (b) the semi-Poincare recurrence Tp/2 ~ 115000. Grid 1200^ 
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VI. APPENDIX : QUANTUM TURBULENCE 
SPECTRA ON 3072^ GRIDS 

Here we consider some detailed spectral simulations on 
3072^ grids, for both winding number n = 1 and n = 2 
straight line vortices. These simulations were run to time 
tmax = 48000. For n = 1 vortices, Ekin{0) / Equ{0) = 
14.0, while for n = 2 vortices ^/cm(0)/^gn(0) = 34.7 
initially. By about t = 10000 this ratio between the 
total kinetic and quantum energies has asymptoted to 
~ 1, with larger fluctuations for n = 1 winding number 
vortices. 

In breaking down the initial kinetic energy into its com- 
pressible and incompressible components. Initially nearly 
all the kinetic energy is incompressible 

n = 1 : ^^°r(0)/^ir™"(0) = 0.0076, (29) 
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n = 2 : i^Hr(0)/-E™''(0) = 0.0052 



(30) 



For winding number n = 1, the compressible energy in- 
creases quite rapidly and becomes equal to the incom- 
pressible kinetic energy around t = 11000, while forn = 2 
vortices ^^^^^ - by t = 4000. By t = 48000, for 

n = 1 E^^^ /E^^^^^^ ~ 2.3 but the state is still evolv- 
ing, while for the more turbulent n = 2 vortices, this ra- 
tio reaches steady state by t = 20000 and then fluctuates 
gently around ^^rV^ir"'' 4.0. The initial spectra 
are plotted in Fig. [l3] for quantum, incompressible and 
compressible energies. 

By t = 10000, the quantum energy spectrum and com- 
pressible kinetic energy spectrum are quite similar. Fig. 
[T4| especially for /c > 10 One is now starting to see 
the development of a triple cascade spectrum for both 
the compressible kinetic energy and the quantum energy, 
while the incompressible energy spectrum has what ap- 
pears to be a backward pulse propagating (in time) to 
smaller k from the large k regions. Moreover, there is 
the indication that around the wave number k = 300, at 
which the compressible and quantum energy spectra en- 
ter their third spectral region, there is a distinctt change 
in the incompressible kinetic energy spectral exponent. 
This is further reinforced by examining the spectra at 
t = 48000, Fig. [15] Note the similarity between the 
spectra at t = 10000 and t = 48000. For winding num- 
ber n = 1, the compressible kinetic energy spectrum is 
quite noisy in the large k > 300 region, and this noise is 
considerably suppressed for quantum turbulence driven 
by winding number n = 2 vortices. 

Using linear regression in fixed /c- windows, we com- 
pute the time averaged energy spectral exponents and 
their average deviation for both the n = 1 and n = 2 
winding number vortices. The time average is performed 
from 60 spectral exponents and their standard deviations 
determined from times t = 24000 to t = 48000 in steps 
of At = 400. 



TABLE II. The time- averaged spectral energy exponents, k 
the various /c-bands for winding number n = 1. Grid 3072^. 



K-BAND 


Energy 


< a> 


± < a > 


15 < /c < 90 


Incomp. K.E 


3.19 


0.17 


180 < /c < 280 


Incomp. K.E 


3.134 


0.016 


350 < /c < 900 


Incomp. K.E 


2.9916 


0.008 


15 < /c < 90 


Quantum 


1.95 


0.23 


180 < A; < 280 


Quantum 


7.810 


0.087 


350 < A: < 900 


Quantum 


3.054 


0.060 


15 < /c < 90 


Total K.E 


1.92 


0.23 


180 < A; < 280 


Total K.E 


8.606 


0.084 


350 < A: < 900 


Total K.E 


3.047 


0.096 



A representative plot is shown in Fig. 16 for the spec- 
tral regions in the total kinetic energy at time t = 48000 
for winding number n = 2 vortices, and for the twin spec- 
tral regions for the incompressible kinetic energy around 
k - 300: 

TABLE III. The time-averaged spectral energy exponents, k~^, 
in the various A-bands for winding number n = 2. Grid 3072^. 



K-BAND 


Energy 


< a> 


± < a > 


15 < /c < 90 


Incomp. K.E 


2.72 


0.21 


180 <k< 280 


Incomp. K.E 


3.377 


0.021 


350 < A: < 900 


Incomp. K.E 


3.006 


0.011 


15 < A; < 90 


Quantum 


1.70 


0.17 


lSO<k< 280 


Quantum 


7.802 


0.058 


350 < A: < 900 


Quantum 


3.033 


0.018 


15 < A; < 90 


Total K.E 


1.66 


0.17 


180 < A: < 280 


Total K.E 


8.527 


0.064 


350 < A: < 900 


Total K.E 


3.042 


0.020 



One can make some conjectures as to physics behinds 
these exponents. One possible conjecture is that for large 
k > 350 one is basically looking at the intravortex spec- 
trum of a 'single' vortex as the incompressible kinetic 
energy spectrum is a very robust k~^. Interestingly, the 
compressible and quantum energy spectra also exhibit a 
marked k~^ for k > 350 - which is not a characteris- 
tic of the spectrum of a single linear quantum vortex. 
Around k ~ 350, both the compressible and quantum 
energy spectra lift sharply away from the k~^ spectrum 
as we move to smaller k. Around this sharp rise, the 
incompressible spectrum moves significantly away from 
k~^ to k~^-^^^ for k < 280 - at least for the winding 
number n = 2 case. (For n = 1, the quantum turbu- 
lence is more subdued and this could account for the 
lower exponent of a = 3.134 in the incompressible ki- 
netic energy). Whether this higher exponent is due to 
quantum Kelvin waves on the vortices needs further in- 
vestigation. The small k spectral region is also of some 
interest. It is tempting to identify this region with the 
classical Kolmogorov k~^^^ energy spectrum - as the to- 
tal kinetic energy spectrum of our compressible quantum 
system exhibits a k~^-^^ spectrum, for the winding num- 
ber n = 2 vortices. It must be remembered that when 
dealing with Fourier transforms we have lost all spatial 
information per se, and so the large k- regions are picking 
up information on physical processes with small spatial 
scale separation. While intravortex physics is one such 
effect, there will be effects from the propagation of wave 
fronts and shocks since the BEC gas is compressible. 
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t=0,(n=l) 



t=0,(n = 2) 




100 



500 1000 




(a) n = 1 



(b) n = 2 



FIG. 13. Initial energy spectra for (a) winding number n = 1 vortices, and (b) winding number n = 2 vortices. Blue circles 
incompressible kinetic energy, red squares - compressible kinetic energy, gold diamonds - quantum energy. Grid 3072^ 



t=10000,(n = 2) 



t= 10000, (n=l) 




k k 

(a) n = 1 (b) n = 2 

FIG. 14. Energy spectra at t = 10000 for (a) winding number n = 1 vortices, and (b) winding number n = 2 vortices. Blue circles - 
incompressible kinetic energy, red squares - compressible kinetic energy, gold diamonds - quantum energy. The secondary peak at /c ~ 40 
for n = 2 vortices in the incompressible kinetic energy spectra appears to be like a backward pulse propagating from much larger k at 
earlier times. Grid 3072^ 
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t=48000,(n = 2) 



t=48000,(n=l) 




k k 

(a) n = 1 (b) n = 2 

FIG. 15. Energy spectra at t = 48000 for (a) winding number n = 1 vortices, and (b) winding number n = 2 vortices. Blue circles - 
incompressible kinetic energy, red squares - compressible kinetic energy, gold diamonds - quantum energy. A triple cascade is quite evident 
in both the quantum and compressible kinetic energy spectra. These two spectra only deviate around the transition from the medium k 
to large k cascade, i.e., around k ~ 300. Grid 3072^ 



t = 48000 Total Kinetic Energy , n = 2 
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(b) n = 2 



FIG. 16. Energy spectra at t = 48000 for winding number n = 2 vortices for (a) the total kinetic energy spectrum Ekin{k), and (b) 
the incompressible kinetic energy spectrum E^jJ^^"^^ (k) . The spectral exponents for the E]^in(k) ^ k~'^ are: small-/c region, a = 1.64, 
for the medium-/c-range o; = 8.10 while for the large-/c region a = 3.06. For (b) the incompressible kinetic energy spectrum, we show the 
regression fit for medium-fc range with a = 3.38 (red line) and large-fc region with a = 3.01 (green line). Grid 3072^ 



